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A 2-D Riemann problem is designed to study the development and dynamics of the slow shocks 
that are thought to form at the boundaries of reconnection exhausts. Simulations are carried out for 
varying ratios of normal magnetic field to the transverse upstream magnetic field (i.e., propagation 
angle with respect to the upstream magnetic field). When the angle is sufficiently oblique, the 
simulations reveal a large firehose-sense (Py > P±) temperature anisotropy in the downstream 
region, accompanied by a transition from a coplanar slow shock to a non-coplanar rotational mode. 
In the downstream region the firehose stability parameter e = 1 — fJ>o(P\\ —P±)/B 2 tends to lock in to 
0.25. This balance arises from the competition between counterstreaming ions, which drives e down, 
and the scattering due to ion inertial scale waves, which are driven unstable by the downstream 
rotational wave. At very oblique propagating angles, 2-D turbulence also develops in the downstream 
region. 

I. INTRODUCTION 

Following the publication of the MHD reconnection scenario of Sweet and Parker [TJ [2] , Petschek [3] noted that 
a pair of back-to-back slow shocks bounding the reconnection outflow could significantly raise the efficiency of the 
process by acting as a transition between the inflowing and reconnected outflowing plasma. Ion heating by these 
Petschek-reconnection-associated slow shocks is one of the mechanisms that has been proposed for solar flares [4] [5] 
and the solar wind. However, although in-situ observations of slow shocks in the solar wind exist [BUTT], they are 
relatively rare, suggesting that the MHD picture of Petschek may not tell the complete story. Unanswered questions 
remain concerning the kinetic structure of such shocks in a collisionless plasma and the associated mechanisms leading 
to particle heating. In previous work, kinetic slow shocks were studied numerically in hybrid codes by initializing 
the system with the slow shock jumps predicted by MHD [12 or, later, by the piston [13 and flow-flow methods 
[13] . Some of the main focuses of these works were the backstreaming beam-driven electromagnetic ion- ion cyclotron 
instability (EMIIC) p~3j[T3], which has been suggested to be the cause of the nonsteady behavior of slow shocks, and 
the formation and damping of downstream large amplitude dispersive wavetrains p~2j [15] . Recently, the dissipation 
due to electrons and beam-excited kinetic Alfven waves was studied by the piston method in particle-in-cell (PIC) 
simulations [T6] [T7]. 

In PIC simulations of reconnection, the plasma downstream of the X-line exhibits large firehose-sense (Pn > Pj_) 
temperature anisotropics, as shown in Fig. [I] (This simulation was discussed previously [18]). Panel (a) shows the 
out-of-plane electron current, with the X-point at (x/di, z/di) ~ (0,20), where d{ is the ion inertial scale. In panels 
(b) and (c), we see turbulence in the B x component correlated with the unstable firehose region; a cut of the firehose 
stability parameter e = 1 — /J>o(P\\ — P±)/B 2 at z = — 35di is shown in (d). The firehose-sense temperature anisotropy 
is notable since in-situ observations of the solar wind clearly show that the proton temperature anisotropy is bounded 
by the marginal firehose and mirror mode stability boundaries [T9] , Hence, it is of interest to more closely study the 
temperature anisotropy distribution across the reconnection exhaust far downstream from the reconnection site, in 
order to understand the effect of the self-generated temperature anisotropy on the propagation and steepening of the 
slow shock. 

The Petschek theory of reconnection [3 predicts that the reconnection exhaust will be bounded by a back-to- 
back pair of standing switch-off slow shocks. However, no clear signature of Petscheck shocks has been seen in PIC 
reconnection simulations, such as in Fig. 0HH], or hybrid simulations [20] [21]. This may be due to the relatively 
small domain sizes in the shock normal direction (e^ in Fig. [T]) , although other reasons have been discussed in the 
context of large-scale hybrid reconnection simulations [20 , 22J. To address this issue, we perform 2-D PIC simulations 
that extend the simulation size in the normal direction to ~ 800^ by ignoring the X-line and instead examining the 
conceptually simpler Riemann problem. In the Riemann formulation, waves propagate away from an interface of two 
different uniform states, such as the two sides of a reconnection symmetry line. This set-up more closely resembles the 
reconnection outflow exhaust than that produced by other methods of generating shocks. Similar Riemann problems 
have been carried out in 1-D resistive-MHD [15], including analyses of asymmetric states and the effect of a guide 
field, and in 1-D hybrid simulations [22 . 2-D Riemann problems have also been carried out in hybrid simulations [23] 
where the 2-D downstream turbulence appeared to diminish the downstream wavetrains associated with switch-off 
slow shocks. Upstream perpendicular heating by EMIIC and the subsequent excitation of Alfven/ion cyclotron waves 
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have also been studied in similar 2-D Riemann problems [24j [25] . 

In Sec. II of this paper we introduce our simulation model and the initial set-up of the Riemann problem. In 
Sec. Ill we discuss the general profiles of a run with Obn — 75° (the angle between the upstream magnetic field and 
the normal direction (e^)). Section IV points out that the counterstreaming ions drive e down (increase the firehose- 
sense temperature anisotropy) in the downstream region. In Sec. V we show that a more oblique angle results in a 
lower e at downstream region. The structure of the magnetic field performs a transition from a coplanar decrease to a 
non-coplanar rotation at e ~ 0.25, which differs from the traditional slow shock transition with dispersive wavetrains. 
In Sec. VI the stability of the downstream rotational wave is studied with numerical experiments. The tendency for 
a spatially modulated rotational wave to radiate c^-scale waves is identified. The resulting c^-scale waves counter- 
balance the e decrease driven by the counterstreaming ions. Finally, we summarize the results and discuss potential 
implications in Sec. VII. 

II. SIMULATION MODELS AND DETAILS 

Our PIC simulations use a narrow computational domain, l z x l x = l.Qdi x 1638.4^ to capture the nonlinear wave 
propagation (mainly the slow shock pair in Petscheck's reconnection model) far downstream from the reconnection 
site. The strategy is to use time as a proxy for space in order to reduce the computational burden. The simulations 
presented here are two-dimensional, i.e., d/dy = 0, and periodic in the z — x plane. The initial equilibrium consists 
of a double Harris-like current sheet (although we only focus on a single current sheet) superimposed on an ambient 
population of uniform density n a : 

B z = B Zja ttmh(x/wi)] n Pi e = n h sech 2 (x/wi) + n a , (1) 

where B z ^ a ,nh,n a are constants, the subscript "a" stands for the asymptotic (far upstream) values, "h" stands for 
Harris and Wi is the initial half- width of the current sheet. We initialize both the Harris plasma and background 
plasma with an isotropic Maxwellian distribution. Unlike the initial set-up for reconnection, we begin with a constant 
normal field B x . Although the initial total pressure is balanced, the existence of B x causes a tension force that drives 
wave propagation away from the current sheet in the ^-direction. The B x ^ a /B z ^ a = 0.1 case {i.e., Obn ~ 83° with 
cosObn = B Xja /B a ) corresponds to a reconnection exhaust with a normalized reconnection rate of 0.1. 

This simulation box can be pictured as existing in the reconnection outflow frame. The distance from the recon- 
nection site is estimated as Ca z x time, where Ca z = Bz,a/ V^o^i^a is the Alfven speed based on the reversed 
component of the field. When l z is small our simulations are essentially 1-D, although we do perform runs with larger 
l z to investigate the possibility of developing 2-D turbulence. Waves (e.g., fast, intermediate, and slow modes in the 
fluid model) will propagate away from the central discontinuity at their characteristic speeds, and can steepen into 
shocks, spread into rarefactions, or maintain their initial shapes based on their own nonlinear it ies. 

In our particle-in-cell code p3d [26] , the electromagnetic fields are defined on gridpoints and advanced in time with 
an explicit trapezoidal- leapfrog method using second-order spatial derivatives. The Lorentz equation of motion for 
each particle is evolved by a Boris algorithm where the velocity v is accelerated by E for half a timestep, rotated by 
B, and accelerated by E for the final half timestep. To ensure that V «E = p/eo a correction electric field is calculated 
by inverting Poisson's equation with a multigrid algorithm. 

The magnetic field is normalized to the asymptotic magnetic field B ai the density to the as ymptotic d ensity n a , 
velocities to the Alfven speed Ca = B a j \±§min a , lengths to the ion inertial length di = y 'mil 'jaon a e 2 , times to 
the inverse ion cyclotron frequency = rrii/B a e, and temperatures to rriiC\. Other important parameters are 
Trii/m e = 25, c/Ca — 15, n a = 1, = 1.5, B a — 1, and the asymptotic value of initial T^ e = 0.1, which implies 
that f3 a = 0.4. The initial electron temperature is uniform, while the ion temperature varies so as to ensure pressure 
balance in the x-direction. We take the time step At = 0.0025 and grid size A = 0.025. We usually take Wi = di 
since the thickness of the dissipation region during reconnection is on the di scale. There are ~4x 10 8 particles in a 
single run. Table I gives further details of the various runs. 

Runs a, b, c, d, f, g, h and k (all have l z = 1.6di, except the Obn = 83° case) will be further discussed in this 
work. Even though the 75° case with a larger l z = 6Adi (Run i) shows downstream 2-D turbulence, the evolution is 
quite similar to the narrow Run f. As can be seen in TABLE I, 2-D turbulence tends to occur for oblique (Obn > 75°) 
cases when l z is large enough. Curves plotted throughout the rest of this paper are quantities averaged in the z- 
direction, since most runs discussed here do not have any significant variation in the z-direction. The 2-D turbulence 
of Obn — 83° (Run k) will also be discussed. 
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III. GENERAL FEATURES OF THE 75° CASE 

A representative case o{Qbn = 75°, Wi = Idi, B g = (Run f ; hereafter referred to as the "75° run") at time 200/ft c i 
is documented in Fig. [2] As soon as the simulation begins, a pair of fast rarefaction waves propagate out from the 
discontinuity with speed ~ I.ICa, while the slow shocks with their downstream rotational waves have speed ~ 0.15CU- 
As in ideal MHD with a symmetric initial condition and zero guide field, a pair of switch-off slow shocks are expected 
to follow the fast rarefaction waves, as shown in Fig. |2|b) [15]. A switch-off slow shock (i.e., the strongest slow shock, 
whose tangential magnetic field vanishes downstream of the shock) propagates at the upstream intermediate speed, 
while the downstream linear slow mode speed equals the downstream linear intermediate speed. The linear slow mode 
and intermediate mode are known to become degenerate at parallel propagation when the plasma (3 (plasma thermal 
pressure/magnetic pressure) exceeds 1. The downstream rotational waves are often identified as dispersive wavetrains. 
The essential physics of this wavetrain can be described by a two-fluid model [27 , or by Hall-MHD. 

Further details of the slow shock pair are shown in Fig. [3] The most significant feature differing from the ideal MHD 
model is the presence of a large downstream temperature anisotropy e = 1 — p,§{P\\ — P±)/B 2 shown in panel (a). The 
corresponding magnetic field structure is shown in panel (b), where the left-hand-polarized rotational wave is clearly 
seen (the polarization will be discussed further in the hodograms of Fig.[5|. The value of e drops from 1.0 upstream of 
the slow shock to ~ 0.25 around the nearly const ant- magnitude rotational waves found downstream. The anisotropy 
factor e affects the strength of the tension force, which in the fluid theory is proportional to e(B • V)B//i - When e 
is positive, the magnetic field has a restoring tension force, while a negative value makes the tension force operate in 
the opposite way, driving the firehose instability. Viewed another way, the phase speed of an intermediate mode is 
Cj = \JsCa cos(0bn)- Therefore, as e drops the intermediate mode becomes slower, or even stops propagating, going 
firehose unstable for e < 0. The x-direction heat flux Q x = j d 3 v(^m5v 2 6v x )f , where 5v = v— (v), is also documented 
in (a). The heat flux peaks inside the transition from upstream of the slow shock to the downstream rotational waves, 
and then becomes negligible. This fact is used in the follow-up paper [28]. In panel (c), the parallel ion temperature 
increases sharply in the weak field region, while the perpendicular ion temperature is nearly constant. The electrons 
are nearly isotropic across the shock. In panel (d), the parallel plasma pressure and perpendicular plasma pressure 
are shown, and the nearly constant normal direction pressure balance P x + B 2 /2po indicates the absence of fast modes 
in the reversal region. Panel (e) shows the associated variations in j3 and the local Obn and panel (f) documents the 
density profiles. In order to pin down important kinetic effects not included in the MHD model (such as the large 
temperature anisotropy), the black dotted curves in each panel show the predicted jumps and positions of slow shocks 
in the ideal MHD version of this global Riemann problem [T5] . 

Within the MHD predicted switch-off slow shock (SSS) jump shown in panel (b), a coplanar transition decreases 
the SSS upstream B z to ~ 0.5. After this, the magnetic field structure rotates in the non-coplanar direction and 
exhibits nearly constant \B\ inside the downstream rotational waves. The amplitude eventually drops to the value of 
B x (i.e., the tangential magnetic field vanishes) in the center, as the symmetry of the initial conditions demands. The 
coplanar transition is recognizable as a slow shock transition, where the major enhancements of the temperatures, 
pressures, densities and decrease in \B\ occur. A similar step-like decrease in \B\ due to a slow shock was also noted in 
the downstream of a large-scale hybrid reconnection simulation [20] . The constancy of the total magnetic field inside 
the downstream rotational wave suggests an intermediate-wave-like behavior. 

IV. THE SOURCE OF TEMPERATURE ANISOTROPY: ALFVENIC COUNTER-STREAMING IONS 



Using the Walen relation Vt,d — Vt,u = ±\/ p u ^ul 'Po(Bt,d/ Pd — &t,u/Pu) fc> r switch-off slow shocks or rotational 
discontinuities ("t" for tangential, "u" for upstream and "d" for downstream), an outflow in the z-direction with 
Alfvenic velocity, B ZjU / ^popu ~ Caz = B ZjCL / ^po m i n a is predicted [29]. The energy source of the downstream 
outflow is the difference in the tangential component of the magnetic field across the discontinuity. If we jump to 
the outflow frame (the de Hoffmann- Teller frame), there will be inflowing Alfvenic streaming ion beams from both 
discontinuities along the downstream magnetic field as observed by Gosling in the solar wind [30 and in kinetic 
reconnection simulations [18, 20] [2TJ (3TJ [32] . These counter-streaming ions cause an enhancement in the downstream 
parallel ion temperature and, therefore, the temperature anisotropy. 

In the phase space of the 75° case at time 200/1^ (see Fig. [4| a signature of the counter-streaming beams is not 
obvious in the downstream region (since the ion distribution does not peak in the upper and lower parts of a single 
wave oscillation), perhaps because of the large-amplitude rotational wave. However, Alfvenic backstreaming ions 
in the z-direction (close to the parallel direction in the upstream region) are observed in Fig. [4] The time-of-flight 
effect (faster ions escape farther upstream) slowly broadens the transition region of slow shocks with time. The 
nearly uniform electron distribution in all directions is due to their high thermal conductivity and much lighter mass 
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compared to the ions. 

V. TEMPERATURE ANISOTROPY VS. PROPAGATION ANGLES 

In order to understand how the temperature anisotropy varies with other parameters, we perform runs with different 
upstream angles Obn- Fig. [5] documents the results of Runs a, b, c, d, f and k. From the first column, the downstream 
s tends to lower values in the more oblique cases and the plasma becomes turbulent once e is comparable to or 
lower than ~ 0.25. We can estimate at which Obn the downstream e will drop below 0.25 as follows. In the cold 
plasma limit, the temperature anisotropy due to the Alfvenic counter-streaming ions at the symmetry line is (where 
P|| - P ± - n a miC\ z and \B\ ~ B XiCL ), 

B 2 

Sdown ~ 0.25 - 1 - -f^ tan 2 B jv, c ~ 0.75 -> BN ,c ~ 40° (2) 

This argument qualitatively shows the tendency to develop stronger firehose-sense temperature anisotropics for higher 
obliquities. The difference between 0bn,c an d the observed value of 60° is probably due to the simplified assumptions, 
such as a cold streaming plasma. 

In the second column of Fig.[5j the corresponding magnetic structures are shown. When the obliquity is large enough, 
especially when e < 0.25, the downstream magnetic field rotates into the out-of-plane direction and becomes turbulent. 
Combined with the hodograms in the third column, we can deduce that the dominant downstream rotational waves 
are all left-handed (LH, counter-clockwise in our hodogram). When the wavelength of the primary LH wave is large, 
as in the Obn = 60° case, its front part breaks into finer right-handed (RH) waves with scale ~ 6d$. In the 75° case, 
the scale of the primary LH wave is already as small as 6<ii, and so it is more stable than the 60° case, albeit still 
turbulent. (A 75° case with a wider initial current layer is discussed in the next section. It exhibits wave-generation 
phenomena similar to the 60° case). In the 83° case, we observe RH small-scale waves in front of the downstream 
primary LH wave. 

For comparison, the dotted curves in the second column of Fig. [5] are the predicted B z structure from MHD theory 
[T5] . The overall predictions agree well for the oblique cases (see, for instance, the upstream B z of the slow shocks in 
the 60° and 75° cases), although the reflected weak fast rarefactions from our boundary have caused a discrepancy 
in the slow shock upstream B z for the 83° case. In less oblique cases, the intermediate and fast characteristic speeds 
approach one another just upstream of the switch-off slow shock according to MHD theory. Therefore there is no clear 
separation between the slow shocks upstream and the fast rarefactions, as can be seen in the simulations. We treat the 
place where the inflow speed V x (not shown) starts to decrease as the upstream of the slow shocks, which corresponds 
to the beginning of the LH rotational wavetrains in the 30° and 45° cases. Their upstream will hence correspond to 
B z ~ 0.15 and B z ~ 0.2 respectively. Therefore the stable small amplitude rotational waves in the Obn = 30° and 
45° cases are more similar to the conventional dispersive stationary downstream wavetrains, which immediately follow 
the slow shock upstream. (We note that the model in Lin and Lee [15 approximates the rarefactions by replacing 
the energy jump condition in the Rankine-Hugoniot relations by [Pp~ 7 ] = 0, arguing that the entropy across a weak 
rarefaction does not change. This is only valid for weak rarefaction waves, but the overall tendency as the propagation 
angle becomes more parallel should be in the correct sense). 

An interesting feature in the oblique cases is the coincidence between the start of the primary LH magnetic rotation 
and the location where the anisotropy parameter e ~ 0.25. The anisotropy parameter seems to be locked to this 
critical value for long periods of time, as is shown in Fig. [6j Moreover, this is also the location where turbulence 
develops. Given these coincidences the following questions naturally arise. If the rotational wave is really a normal 
dispersive wavetrain, why does it appear in the middle of the MHD predicted switch-off slow shocks? Why do the 
dispersive waves not start directly from the slow shock upstream? Is it not more similar to a new transition at 
B z ~ 0.3 for B N = 60°, or B z ~ 0.5 for B n = 75° and 83°? What is special about e = 0.25? The importance of 
e 0.25 is shown for different obliquities in Fig. [7[ Are there other instabilities associated with this anisotropy value? 
Or is it due to the nonlinear structure of a system with a large temperature anisotropy that cannot be explained by 
ideal MHD? 

We describe a possible theoretical explanation for why e ~ 0.25 in a follow-up paper [28 . In short, e = 0.25 
represents a transition where the slow and intermediate mode speeds become degenerate. Unlike the conventional 
picture, the downstream rotational waves can not then be explained by slow dispersive waves, but instead take the form 
of rotational intermediate modes. The coplanar decreasing part of the magnetic field together with the non-coplanar 
rotational part will later be identified as a single nonlinear wave, called a compound SS/RD wave [28] . 
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VI. THE DOWNSTREAM TURBULENT WAVES AND PARTICLE SCATTERING 



The region downstream of oblique (6bn > 60°) slow shocks is turbulent. High wavenumber waves are continually 
excited whenever the B z component begins to rotate into the out-of-plane direction, as can be seen for the 60°, 75° 
and 83° cases in Fig. [5| and more clearly in the evolution of the B field for the case 75°, W{ = 10di (Run g) in 
Fig. [8] The downstream LH rotational parent waves break into A^ ~ 6<ii-scale waves. The large oscillation in e at 
later times is due to the small magnetic field magnitude near the symmetry line, where both B z and B y vanish. The 
particle scattering associated with these small-scale waves plays a role in counter-balancing the decrease in e due to 
the streaming ions and keeping the temperature anisotropy around the value 0.25, as is seen in the time evolution of 
the oblique cases in Fig. [6| 

In order to understand the downstream turbulent <^-scale waves, we tried to systematically pin down the possible 
driver and energy source via numerical experiments. We separately checked potential energy sources in the shock 
simulations, including temperature anisotropy, counterstreaming beams, rotational parent waves (such as the larger 
LH wave in Fig. [8]) to determine which factors are responsible for generating the <^-scale waves. 

We carried out spatially homogeneous simulations with an initial wave structure of the following form: 

B z = B cir cos(2ttx/A p ) + ^oblique; B y = P x B cir sin(27rx/A p ) (3) 

When .B c i r ^ 0, there is a rotational field with P = +1 for LH, —1 for RH, and for planar polarizations. The 
constant oblique controls the obliquity of this circularly polarized wave and provides a spatial modulation in the 
magnitude of the total magnetic field. The general expression for the initial ion distribution is 

f oc exp ( mi[V]] - U)2 mv2 A I expf m ^H +M ) 2 m * v2 A U) 
U H 2T.ii 7Uy! +eX H 2Ii|| T i± ) W 

which has bi-Maxwellian counterstreaming beams for u ^ 0. Note that Ti||,eff = T^i + m^ 2 , so both u ^ and 
Ti\\ 7^ Ti± can contribute to the temperature anisotropy e. The initial plasma density varies so as to ensure a constant 
value of P± + B 2 /(2fio). Since the small-scale waves that interest us do not induce variation in the z-direction, they 
are intrinsically 1-D waves along the x-direction. The common parameters used here are a domain size of 1.6 x 51.2, 
with grid 64 x 2048, X p = 51.2^, uniform T^ e j_ = T e u = 0.15, B x = 0.25 and n^ e = 1.5 at x = 0. These parameters 
are meant to represent those observed in the downstream of the 75° case. 

Run 1 of Table II is a representative example of the downstream structure seen in the Riemann simulations. The 
initial obliquely propagating LH polarized waves with Alfvenic counter-streaming ions along the local magnetic field 
break into ~ 6di small-scale waves, as shown in Fig. |9ja). In Run 2, we replace the counter-streaming ions with a 
bi-Maxwellian plasma with the same effective parallel temperature, therefore the same temperature anisotropy e, and 
find similar wave-generation phenomena. In Runs 3 and 4, we remove the obliquity (B z oblique) from Runs 1 and 2 
respectively. No <^-scale waves are excited within a time of 100/r2 C i in Run 3, as is shown in Fig.[9^b). This indicates 
the importance of spatial modulation for the development of the turbulence. In Run 5, we replace the beams or 
anisotropic plasma of Runs 1 and 2, respectively, with a Maxwellian isotropic plasma. Waves with scale ~ 6di are 
excited and are shown in Fig. [9jc) . However, if we further remove the initial out-of-plane magnetic field of the parent 
wave (Run 6), no small-scale waves appear; see Fig. |9jd). This indicates that a necessary condition for producing 
these small-scale waves is the existence of circularly polarized parent waves. From this suite of runs we conclude that 
the presence of a spatially modulated rotational wave is the major driver of descale coherent waves seen downstream 
of the shock (see Fig. |8jb)). Although not shown, we see similar behvior for RH (P = —1) parent waves. 

A non-modulated, constant-amplitude, circularly polarized Alfven wave is a known stationary solution of the MHD 
equations. Although it is subject to a long- wavelength modulational instability under some conditions [33], simulations 
show that it is also stable in a collisionless plasma. Spangler [34] explored a similar wave-generation phenomena, as 
well as the formation of solitary waves with a system with initial conditions similar to our Run 5. He cast his 
system in terms of the Derivat ive- Nonlinear- Schrodinger-Equat ions, which are simplified MHD equations that possess 
the MHD nonlinearity and dispersion terms. This suggests the importance of both the nonlinearity and the Hall 
dispersion term for these d^-scale waves. We also note that the time of onset of these waves is proportional to X p and 
inversely proportional to the amplitude B c [ r , which are closely related to the steepening time scale of finite amplitude 
waves, T s ~ (r wave /27r)(B / AB), with r wave oc X p (the gradient scale of the finite amplitude wave) [35]. Although a 
firehose-sense temperature anisotropy (e < 0) would weaken the nonlinearity, the spatial variation of e enhances the 
nonlinearity, as will be shown in the follow-up paper [28] . 

Runs 3 and 4 suggest that temperature anisotropy and Alfvenic beams are not, by themselves, sufficient for 
generating the small-scale waves seen in the downstream region. The time scale of the wave- generation is faster than 
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beam/anisotropy driven waves from the electromagnetic ion- ion cyclotron instability (EMIIC), such as Alfven/ion 
cyclotron (AIC) waves [36 , or kinetic- Alfven waves (K AW) [37] . As is the case for the growth rate of the firehose 
instability 7 2 ~ (— sk 2 C\)/(l + C\/c 2 ) [38], the growth rate is small when Ca (based on the ambient B field) is 
small. However, these beam/anisotropy-driven waves are still potential players in the downstream turbulence at late 
time. By comparing Fig. [9|a) to (c), we see that the interplay of the steepening process and the free-streaming beams 
makes the wave more turbulent. In Fig. ~~ ~~ 
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to Fig. [9ja)) decreases. This confirms the a 



the ion temperature anisotropy — Ti± of Run 1 (which corresponds 
ibility of these smaller-scale dispersive waves to scatter ions. 
The polarization of a linear wave can be determined by the phase between the tangential magnetic field variations 
SB Z and SB y . In the fluid model iSB z /6B y = (C 2 — (u / 'k) 2 ) / \uodiC Ucos(6 ] bn)) [H.- A wave has a LH polarization 
when (uo/k) 2 > Cf , and RH otherwise. As observed in our oblique shock simulations (discussed in Sec. the 
small-scale waves at the upstream of the rotational front are mostly RH, while the downstream waves are mostly LH. 
This suggests that the primary rotational front propagates at the local intermediate speed (as also measured in the 
simulations). We do not address which modes are responsible for these fine-scale waves here, since the polarization of 
the linear mode in kinetic theory is very complicated. Both temperature anisotropy and high plasma (3 play roles in 
changing the linear wave properties [39.. With a larger l z and very oblique propagation angle, as in the 83° case, 2-D 
turbulence with X z ~ 6di is excited around the firehose-unstable region. Its signature is most clear in B x , as can be 
seen in the bottom plot of Fig. [6] A similar mechanism is postulated to excite 2-D turbulence in the firehose-unstable 
region seen in the reconnection simulation of Fig. [TJc) (and, perhaps, the 2-D turbulence previously reported [23]). 
This 2-D turbulence could be driven by temperature anisotropy (firehose-like) or by the sharp front of the primary 
rotational wave. 



VII. SUMMARY AND DISCUSSION 



We have studied the temperature anisotropy distribution across slow shocks with different obliquities in PIC sim- 
ulations. An abnormal transition and an anisotropy e = 0.25 locking phenomena downstream of the MHD predicted 
switch-off slow shocks is documented. The Alfvenic counter-streaming ions serve as the driver for decreasing e (in- 
creasing the firehose-sense temperature anisotropy) in the downstream region, while downstream ^-scale turbulent 
waves scatter particles and raise e. This dynamical balance makes the downstream e plateau at a value of e = 0.25 
and not the marginal firehose criterion e = 0. The theoretical significance of e = 0.25 will be addressed in another 
work [28 . By means of PIC numerical experiments we show that the turbulent c^-scale waves are radiated from a 
spatially modulated rotational parent wave . 

The Riemann problem for the Obn — 83° case is closely related to that of reconnection exhausts with normalized 
reconnection rates of 0.1. The very center is a firehose-unstable region where a B x variation is observed, as in the 
reconnection simulations. Although we cannot confidently identify the 0.25 plateau in present PIC reconnection 
simulations (since the spatial extension in the normal direction (e x ), is ~ 10^; see Fig. [TJd)), we expect to see the 
signature of a e = 0.25 plateau outside the firehose unstable region in very large kinetic anti-parallel reconnection 
simulations and in-situ satellite observations of anti-parallel magnetic reconnection outflows. 

Compared to slow shocks, fast shocks have been intensively studied (see, for instance, a review article [40]). The 
formation of shocklet and short large-amplitude magnetic structures (SLAMs) in front of earth's quasi-parallel bow 
shock has been observed [41] and simulated [42j[43]. In our slow shock simulations, the development of downstream 
rotational waves produces SLAMs-like structures as in Fig. |8Fb). Since the transition of the magnetic field across 
the front in a fast shock is opposite to that in a slow shock, if there are SLAMs associated with a slow shock, they 
are expected to exist in the downstream region, especially since the downstream of a switch-off slow shock is locally 
quasi-parallel. 

Finally, even though we have not seen super-thermal particles other than Alfvenic streaming ions in our shock 
simulations, it is still an open question whether particles could be accelerated by slow shock associated reconnection 
exhausts via mechanisms recognized in fast shocks such as the lst-order Fermi mechanism [44], diffusive shock ac- 
celeration (DSA) [45, 46 , and those possible injection mechanisms for DSA such as shock-drift [47], or shock surfing 

EH]. 
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FIG. 1: The exhaust from steady reconnection in a PIC simulation. Panel (a): The out-of-plane electron current J ey ] Panel 
(b): B 2 s, where positive values have been set to 0. The colored region is firehose unstable; Panel (c): The magnitude of B x 
showing the turbulence associated with the firehose instability; Panel (d): A cut of e at x/di — —35 (the vertical line in (b)). 
The horizontal lines demark e — 0.25 and e — 0. 




FIG. 2: The evolution of a system with Obn — 75° (Run f). Panel (a): The evolution of B from time — 200/f^. A pair 
of fast rarefactions (FR) propagate out from the symmetry line, followed by a pair of slow shocks (SS). Each curve has been 
shifted so that it intersects the vertical axis at the given time. The time between the yellow curves is 100/n c i; Panel (b): The 
predicted FR and switch-off slow shock (SSS) from ideal MHD theory; Panel (c): The same as (a) but with the vertical axis 
measuring B; Panel (d): The evolution of B z from time — 200 /Q C i- 
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(a) (b) 




FIG. 3: Parameters from the run with Obn — 75° (Run f) at time 200/n c i. Panel (a): Temperature anisotropy e and 
x-direction heat flux Q x ; Panel (b): Magnetic field components; Panel (c): Parallel and perpendicular temperatures (the off- 
diagonal components Ti xy ,Ti XZ ,Ti yz are plotted together in green, denoted as T Q ff, and are small) ; Panel (d): Total plasma 
pressure components and P x + B 2 /2/j j o. Panel (e): The plasma f3 and local Obn — cos -1 (B x /B); Panel (f): Plasma density. 
The dotted curves in each panel are the predicted magnitude and position of the switch-off slow shocks (SSS) from isotropic 
MHD for B z in (b), T in (c), P in (d), f3 in (e), and n in (f). 
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FIG. 4: The phase space of the run with Obn = 75° (Run f) at time 200 /Q C i- From top to bottom the left column shows the 
ion distribution in: V z — x space, where the backstreaming ions from the discontinuities are clearly seen; V y — x space; V x — x 
space. The right column is the electron distribution in V z — x space, V y — x space and V x — x space. 
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FIG. 5: From top to bottom are runs with B n = 30° (Run a) at 100/Q ci , 45° (Run b) at 200/Q ci , 52° (Run c) at 100/Q ci , 60° 
(Run d) at 250/Q C i, 75° (Run f) at 400 /Q C i, and 83° (Run k) at 700 /Q c i- The first column shows the temperature anisotropy, 
and the second column the magnetic field components as a function of x. The third column displays hodograms taken from 
the right half of the simulation domains. The dotted curves in the second column are the predicted magnitudes and positions 
of switch-off slow shocks (SSS) and fast rarefactions (FR) from isotropic MHD theory. 
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FIG. 6: Evolution of e for the case with Obn — 60° (Run d) for equally spaced times between 100 — 500/n c i from lighter grey 
to darker grey in (a), the Obn = 75° case (Run f) for time 100 — 500 /Q C i in (b), and the Obn = 83° case (Run k) for time 
100 — 700/Q C 2 in (c). The bottom is a plot of B x for the Obn — 83° case at time 700/f2 C i showing the 2-D turbulence that 
develops. 
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FIG. 7: The e distributions of runs 6 BN = 60° (Run d) at 500/tt ci , 75° (Run f) at 200/n ci , 83° (Run k) at 700/n ci . (The 
83° case is shifted to the right by 204.8di) 
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FIG. 8: Panel (a): The evolution of e, B 2 and B y for equally spaced times between 50 — 500/f2 c ; (from left to right) in the 
Qbn — 75°, Wi = 10di case (Run g). The downstream larger-scale rotational wave breaks into waves of wavelength ~ 6di. 
Panel (b): A blowup of the downstream B y at time 450 /Q c i- 




FIG. 9: The evolution of B y , B z and B for equally spaced times between — 100/^ C i- The red curve indicates the time 100/^ C i- 
Panel (a): Run 1 with both initial streaming ions and modulated rotational parent wave. Panel (b): The same as panel (a) 
without the initial spatial modulation (Run 3). Panel (c): The same as panel (a) without initial beams (Run 5). Panel (d): 
The same as panel (c) without initial polarization (Run 6). 
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FIG. 10: The evolution of T^y — Ti± for equally spaced times between — 100/n c i(f rom lighter grey to darker grey) of Run 1 
(Fig. (9ja)). The temperature anisotropy of the ions is reduced, which indicates particle scattering is taking place. 



Tables 



TABLE I: 



Run 


Obn 


Wi B g a 


Domain Size (l z x l x ) 


Gridpoints - 


V 2-D turbulence 


a 


30° 


1 


1.6 x 1638.4 


64 x 65536 


X 


b 


45° 


1 


1.6 x 1638.4 


64 x 65536 


X 


c 


52° 


1 


1.6 x 1638.4 


64 x 65536 


X 


d 


60° 


1 


1.6 x 1638.4 


64 x 65536 


X 


e 


60° 


1 


6.4 x 819.2 


256 x 32768 


X 


f 


75° 


1 


1.6 x 1638.4 


64 x 65536 


X 


g 


75° 


10 


1.6 x 1638.4 


64 x 65536 


X 


h 


75° 


1 0.2 


1.6 x 1638.4 


64 x 65536 


X 


i 


75° 


1 


6.4 x 1638.4 


256 x 65536 


V 


j 


83° 


1 


1.6 x 819.2 


64 x 32768 


X 


k 


83° 


1 


6.4 x 819.2 


256 x 32768 


V 



B g is an initial uniform guide field in the y-direction 
'"—>." means "resulting in". 



TABLE II: 


Run B c i r 


Polarization(P) B z , oh u que 




Beams (u) 


. a rri 

1 M\,eff 


— ► € 


-» 1-D di-scale waves (t < 100/Q C i) 


1 0.25 


+ 1 


0.25 


0.15 


0.5 


0.4 


-0.2 - -7 


V 


2 0.25 


+ 1 


0.25 


0.4 





0.4 


-0.2 - -7 


V 


3 0.25 


+ 1 





0.15 


0.5 


0.4 


-2 


X 


4 0.25 


+ 1 





0.4 





0.4 


-2 


X 


5 0.25 


+ 1 


0.25 


0.15 





0.15 


1 


V 


6 0.25 





0.25 


0.15 





0.15 


1 


X 



a "— >" means "resulting in" 



